Non-adiabatic Coulomb effects in strong field ionisation in 
circularly polarised laser fields I: Ionisation rates 

Jivesh Kaushal and Olga Smirnova 

Max Born Institute, Max Born Strasse 2a, 124.89 Berlin, Germany 

Abstract 

We develop the recently proposed analytical R-matrix (ARM) method to encompass strong field 
ionisation by circularly polarised fields, for atoms with arbitrary binding potentials. Through 
ARM, the effect of the potential can now be included consistently both during and after ionisation, 
providing a complete picture for the effects of the long-range potential. We find that the Coulomb 
effects modify the ionisation dynamics in several ways, including modification of (i) the ionisation 
(exit) times, (ii) the initial conditions for the electron continuum dynamics, (iii) the "tunnelling 
angle" , at which the electron "enters" the barrier, and (iv) the electron drift momentum. We derive 
analytical expressions for the Coulomb-corrected ionisation times, initial velocities, momentum 
shifts and ionisation rates in circularly polarised fields, for arbitrary angular momentum of the 
initial state. We also analyse how the non-adiabatic Coulomb effects modify (i) the calibration of 
the attoclock in the angular streaking method, and (ii) the ratio of ionisation rates from p~ and 
p+ orbitals. 

PACS numbers: 32.80.Rm, 42.50.Hz, 33.80.Wz 



I. INTRODUCTION 



Single and double ionisation in circularly polarised strong laser fields is a sensitive probe of 
attosecond dynamics [TH6] . Strong field ionisation is often viewed as electron tunnelling from 
atoms and molecules through the barrier created by the laser field and the core potential. 
The adiabatic approximation, frequently used to describe tunnelling, implies quasi-static 
electric field and zero electron velocity immediately after ionisation (at the tunnel exit). 
This adiabatic picture is used for the interpretation of current experiments in circularly 
polarised laser fields [IHS] within the two-step model. This model merges quantum and 
classical approaches by combining (i) the adiabatic approximation for the quantum ionisation 
step with (ii) the classical trajectories calculation after tunnelling. In this second step, 
an ensemble of classical trajectories is launched outside the barrier; the distribution of 
initial velocities parallel and perpendicular to the direction of the instantaneous laser field 
is centered around zero, as predicted by the adiabatic tunnelling theory. 

However, strictly speaking, in circularly polarised laser fields the tunnelling barrier is 
rotating. This rotation manifests itself in the non-adiabatic electron response, which becomes 
significant in the regime of the Keldysh parameter 7 > 1. Non-adiabatic effects change 
tunelling from essentially one-dimensional, characteristic of the static limit IH E]) to two- 
dimensional. As shown in [8J, for short-range potentials substantial deviations from the 
adiabatic approximation arise already for 7^ ~ 0.5, which also questions the validity of this 
approximation for the long-range core potentials under similar conditions. We note that 
7^ ~ 0.5 is a typical regime for recent experiments with laser radiation around 1600-1300 
nm and systems with ionisation potential Ip ~ 10 eV (see e.g. [9]). 

Here we provide rigorous analytical framework for treating the effects of long-range po- 
tential and laser field on equal footing and include non-adiabatic effects due to the long 
range potential. In particular, we show that non-adiabatic Coulomb effects lead to non-zero 
initial velocity both parallel and perpendicular to the direction of the instantaneous laser 
field even for the central electron trajectory (for which the ionisation rate maximises), when 
it emerges from the classically forbidden region. 

We find that the non-adiabatic Coulomb effects modify the ionisation dynamics in several 
ways, including modification of (i) the ionisation (exit) times, (ii) the initial conditions for 
electron continuum dynamics, (iii) the "tunnelling angle", at which the electron "enters" 
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the barrier, and (iv) the electron drift momentum. We derive analytical expressions for 
the ionisation times, initial velocities, momentum shifts and ionisation rates in circularly 
polarised fields for arbitrary angular momentum of the initial state. We also analyse how 
the non-adiabatic Coulomb effects modify (i) the calibration of the attoclock in the angular 
streaking method, and (ii) the ratio of ionisation rates from p~ and orbitals obtained for 
short range potentials in f8]. 

Our tool is the gauge-invariant, analytical R-matrix-type method (ARM), which we have 
recently developed [lU] for linearly polarised fields. The strength of ARM is the ability 
to treat consistently the effects of long-range potential and the laser field (TOj as well as 
multielectron effects [TT]. The main idea of the method was adopted from the study of 
collision processes and nuclear resonance reactions [12], where the primary purpose was 
to isolate the strongly interacting kernel from the region where these interactions were 
significantly weaker, and can be considered in asymptotic approximation. In this sense, 
the R-matrix approach developed in collision physics meets a crucial requirement of strong 
field physics, in that it can be used to separate the region of the configuration space where 
the Coulomb forces are much stronger than the laser field, from the outer region where the 
Coulomb potential quickly becomes almost negligible compared to the driving laser field. 
We note that fully numerical time-dependent R-matrix approach for strong-field dynamics 
has been developed and successfully applied in [13]. 

In [in], detailed analysis and benchmarking of the ARM method was provided for strong 
linearly polarised laser fields, including the derivation of analytical results for the instan- 
taneous ionisation amplitudes and the sub-cycle ionisation rates for single active electron 
systems. With suitable approximations, the results for the cycle-averaged ionisation rates 
from PU] agree with those obtained by Perelomov, Popov and Terentev (PPT) [13]. For 
hydrogen-like atoms and ions, the PPT rates (with the new correction factor derived in |15j ) 
were shown to be accurate for arbitrary values of the Keldysh parameter. 

The difference between our problem in circularly polarised fields and the case of the 
linearly polarised fields analysed in [TU] lies in the fundamentally two-dimensional character 
of tunnelling, i.e. the time-varying angles between the position vectors, the electron velocity 
Vp(t), and the laser vector potential A(t). As a consequence, certain approximations in [TO] 
that were helpful in deriving physically transparent solution for the linearly polarised field 
are not always adequate in the case of circularly polarised fields. Here, we refine the ARM 
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method and extend it to the strong circularly polarised fields. We show that an appropriate 
choice of the the boundary between inner and outer region allows one to build the hierarchy 
of interactions and show how the long-range effects can be included consistently within the 
iterative approach. 

Our strategy can be summarised as follows. Following pLOj , we introduce the "R- matrix" 
sphere of radius a, which splits the configuration space into the inner and outer regions. 

A. The inner region: In the inner region the Coulomb field dominates, and the effects 

of the laser field on the inner region wave-function can be included in the quasistatic 
approximation. Further approximation, such as using the field-free wave-function is 
justified for fields significantly smaller than k'^, where k = ■^2Jp, Ip is the ionisation 
potential. 

B. The boundary: The Bloch operator is used to "pass" the information about the elec- 

tron wave-function from the inner region to the outer region. The outer region Green's 
function is used to propagate the outer region wave-function from this boundary to 
the detector. The boundary value is given by the inner region wave-function at the 
surface of the "R-matrix" sphere. Boundary matching ensures that the final result 
does not depend on the choice of the boundary value a. 

C. The outer region: In the outer region, the Coulomb potential is weak and can be in- 

cluded in the Eikonal-Volkov approximation (EVA) [16j. It has been shown previously 
|17] . that the Eikonal-Volkov approximation is adequate for the soft-core potentials, 
i.e. for the long-range potentials outside the singularity region. 

D. Propagation in the outer region: involves integration over the surface of the R- 

matrix sphere (6*', cj)') and over all times {t') of "transition" through the boundary. 
Due to the large action S of the electron in the strong laser field, the integrals are 
taken from highly oscillating function P{9' , (f)', t')e*'^^^'''^''*'-' and are accumulated in the 
vicinity of their respective stationary (saddle) points 9'^, 0'^, t'^, defined by the solutions 
of equations Sq/ = 0, S(f,> = 0, and Sf = respectively, where the subscript denotes 
derivative w.r.t. that variable. The action is given by = 5*^^^ + Gc, where the 
SFA (Strong Field Approximation) action is associated with the dynamics in the laser 
field and short range potential, and Gc is the action associated with the interaction 
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with the long-range potential of the core under EVA [16] and describes Coulomb-laser 
coupling [IH]. Since only a vicinity of saddle points contributes to the integral, we 
make Taylor expansion of Gc around SFA saddle points 6's'^\ and After this 
expansion the integral over the sphere surface is calculated exactly. The integral over 
t' is evaluated using the saddle point method. The actual, full saddle point t'^ (shifted 
from due to long-range effects) is found within the iterative approach. Formally, 
non-adiabatic effects in ionisation rates arise due to the deviations from the stationary 
trajectory included via Taylor expansion of Gc- Non-adiabatic coulomb effects also 
manifest itself in the photoelectron spectra and will be considered in our subsequent 
paper [T9] . 



E. Iterative approach to saddle point equation for t': By construction, in the outer 
region Gc presents a perturbation to the SFA action S^^^ and therefore can only 
slightly shift the SFA saddle point t'j^\ which correspond to the stationary SFA action: 
dt'S^^^ = 0. Thus, as a first correction to the saddle point, due to the interaction with 
the long-range potential tT^ = tf^ + AtT\ where Atf^ can be found by iterations 
with respect to Gc- The first iteration includes only linear terms in At'j^^ ~ Gc- In 
our approach we keep only the first-order correction terms consistently throughout. 
The saddle point equation St' = can be expanded around tf^: dt'S^^"^ (t'P) + 

Atf^d^,S^-^ (tf + d,Gc (tf^) = 0, yielding Atf ^ = Note that 

since the SFA action is stationary dfS^^^ ^t'i'^''^ = 0, the shift due to At'j^^ will 
only change the value of the SFA action in the second order wrt Gc- S (^'s^^^ = 
SSFA (^^;{o) j ^ 1^^/(0)^ 0f \ tf + O (Atf . However, tf ^ will contribute to the 
pre-exponential factor P[9\(j)',t') in the integral. 

Below we detail our method and show how it can be used to obtain ionisation amplitudes 
and ionisation rates using both the time-domain and the frequency domain approaches. The 
time-domain approach is technically simpler and allows one to consider temporal dynamics of 
ionisation, including time evolution of electron momentum distributions [19] and ionisation 
rates. 

The paper is organized as follows. Section [TT] introduces basic equations. Section HI 



develops the time-domain approach. Section IIV] discusses the physical picture arising from 



the theory developed in sections [IT]- III In section IV, we describe modifications of the 



ionisation dynamics due to the Coulomb effects. These include (i) Coulomb corrections to 
ionisation times, (ii) initial conditions for electron continuum dynamics and calibration of 
the attoclock in the angular streaking method, (iii) Coulomb corrections to the "tunnelling 
angle", including the Coulomb corrections to the ratio of ionisation rates from p~ and 
p"*" orbitals obtained for the short range potentials in [8j . Section |v] concludes the work. 
Appendix |X] presents additional calculations related to the boundary matching. Appendix [B] 
develops the frequency-domain approach, pioneered in the PPT work for the short-range 
potentials. This approach requires more involved algebra but allows most straightforward 



connection to the PPT results. Appendix [C] extends the time domain method in section [ITT 
to include subcycle dynamics. Appendix D -[E] present miscellaneous calculations. 

II. BASIC EQUATIONS 

Following [To], we introduce the Bloch operator (a) to split the configuration space 
into the inner and outer regions. Parameter a represents the radius of the R-matrix sphere, 
the inner region is inside the sphere, the outer region is outside of the sphere. The standard 
Hamiltonian H including both Coulomb Vc{r) and laser-field interaction VL(t): 

^ = y + Voir) + VLit) (1) 



used in the Schordinger equation: 

.dijj (r, t) 



HiPiv.t) (2) 



dt 

ij{r,t = to)=ijg{r) (3) 

can be modified to 

^d^^ = Hsi^ (r,t) - L±(a)^ (r,t) (4) 

where Hb = H + L^{a). Following arguments developed in [10], we can express the solution 
in the outer region via the solution in the inner region as: 

IV^outW) = I f dt'UB{t,t')L-{a)\^Ut')) (5) 

J to 

where for the outgoing solution, we use L~{a) and the governing equation for the evolution 
operator (Jsitji') is: 

t^UB{t,t') = HB{t)UB{t,t'). (6) 



In our time-domain approach detailed in the next section, we start the analysis from the 
expression for the ionisation amplitude ap(T) = {p\ipout{T)) (see for discussion): 

ap(T) = ^l^dt' J dv' j dv" (p|f/B(T,t')|r')(r'|L-(a)|r")(r'>in(t')) (7) 
Taking into account the explicit form of Bloch operator in coordinate representation: 



d 1 1 , / 



(r'|L-(a)|r") = 5{r - a)5(r' - y")B, 
H (r,t) = 
we can rewrite Eq. ([T]) as follows: 

ap{T)=ij^ dt' j dr' GB{p,T;r',t')6{r' -a)B{a,9',(f)',t') 



(8) 
(9) 

(10) 



where function B{a, 6', 0', t') represents the inner region wavefunction at the partition surface 
r' = a for all times to < t' <T: 



Bia,e,<P,t')= (^|- + i^V^i„(r,t') 



(11) 



and Gb(p, T; r', t') = {p\UB{T,t')\r') is the Green's function for modified Hamiltonian Hb 
for propagating from the boundary r' = a instead of the origin. As shown in [10], the error 
incurred in approximating this exact Green function with the Eikonal-Volkov approximated 
Green's function G^^^(p,T; r',t') defined on the EVA states [16|: 



G'^v^(p,T;r',0 



1 g-»vp(t')-r'-f /Jdri,2(r)gi/^'l/(rz,(r;r',p,t'))g-iGop{ri(T;r',p,t')) 



(27r)3/2 

is exponentially small. In the above expression we have defined: 



rL(r;r,p,t) =r + ^ rfCvp(C) 



(12) 



(13) 



the characteristic trajectory along which the Coulomb correction is calculated as a pertur- 
bation to the Volkov electron [15] and Vp()f:) = p + A(t) is the kinetic momentum. 

III. THE TIME-DOMAIN APPROACH 



We use Eq. (10) for the ionisation amplitude and Eq. (12) to obtain 
a (T) = ^^ r dt' I rfr'e-^^p(*'>'-'e*^''^(P'^''-'*\^^^(P'^''-''*')(5(r'-a)5(a,^',0',t') (14) 



with the Coulomb phase term defined as: 



Gc{r',t',T)= [ rfrf/(ri(r;r',p,t')) (15) 
Jt 

Since the time T of observation is sufficiently far, so that we can consider T — t- oo for all 
practical purposes, we have made the approximation as in [TU], of ignoring the distortions 



of the phase-front from plane wave: Gop [16j in Eq. (12). The SFA phase is given by 



S'^\p,T;r',t') = -\[dTvl{r) + ^t' (16) 
A. Transition through the boundary r' = a 

Function B{a, 6', 0', t') reflects the value of the inner-region wave-function at the boundary 
r' = a. In the inner region the Coulomb field dominates and the effects of the laser field 
on the inner region wave-function ?/'in(r',t') can be included in quasistatic approximation. 
Following [To], the boundary is placed in the asymptotic region Eoa/Ip ^ 1 ^ /ta of the 
ground state wave-function, where k. = Eq is the amplitude of the laser field. In this 

region the difference between polarised and initial wave-function is of the order of ~ Eqo? / k 
|2Uj . Thus, for sufficiently weak fields such that Eq/h? <^ I/o^k^ the inner region wave- 
function ipi^ir' ,t') can be substituted by the field-free bound state wave-function, without 
affecting the boundary matching. This approximation was first used in the PPT method 
|14j . The asymptotic radial part of this wave- function is given by: 

p-nr' 

^.i{r') = CW — {^^r'f^ (17) 
Kr 

Due to the invariance of the boundary term under addition of a function 6o/r', we can choose 
6o appropriately to get: 



5(a,^',0',t') = ( ^ - ) Vni{r'" 



(18) 



for = Q / K. Using Eq. (18) and evaluating the Delta function over r': 

ap(r) = -^ rdt' I rfl]'e-^p(*>--t/*'^-pWe^G^(P'^^-'*V.K«)A^^rnPr(cos0') 

(2vr)3/2 J ^^g^ 
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where N^m = ^ = a(sin6''cos0'x + sin 6*' sin 0' y + cos6''z). We use A{t) = 

—Aq (cos utic + sinuty). 

Note that in the outer region, the long-range interaction of the electron with the core is 
described by the phase term Gc{p,T]r' ,t'), and involves integration of Coulomb potential 
along the electron trajectory in the laser field. The trajectory originates from the point a 



on the boundary at the time t' (Eq. 13) 



B. Iterative approach to solution of saddle point equations 

The calculation of the ionisation amplitude ap{T) involves integration over all starting 
points of the trajectory on the sphere and all the times t' of "transition" through the bound- 



ary r' = a. The integrand of Eq. (19) can be written in the form P{9', 0', t')e'^^^^'''^'''^'\ where 
the prefactor P{6',(j)',t') reflects the value of the inner-region wave-function at the bound- 
ary r' = a, and exponent is given by the electron action. The action S{6',(f)',t') consists 



of two parts: 5* = 5*^^^ + Gq, where 5*^^^ is the SPA action (Eq. 16), associated with the 



ionisation dynamics in the short range potential, and Gc (Eq. 15) is the term responsible 
for the long-range interaction with the core, describing the Coulomb-laser coupling |18j . 

Due to the large value of the action 5* for the electron in a strong laser field, the integrals 
are accumulated in the vicinity of stationary (saddle) points 9'^, 0'^ and t'^, which satisfy the 
equations: 

dS dS dS dvJt') , , 

where the saddle point in time has an additional term afp(t'), which, as we will show in 



Appendix B 1 is a result of propagation from a finite boundary, and comes from the exact 
evaluation of the surface integral. 

To solve these equations, we recall that by construction, in the outer region Gc presents 
a perturbation to S"^^^ and therefore can only slightly shift the SEA saddle points 6*^°^ 0s 
satisfying the equations: 

dQ' ' 90' ' dt' dt' ^ ^ 

Thus, the saddle points for the total action S = S^^^ + Gc can be written as: 

0: = ef) + AC 0l = 0f) + A0'„ C = t'i°) + At; (22) 
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where A^^, A0g, At'^ are the small corrections to the SFA saddle points and can be found 
perturbatively. Subscript "a" in At'^ indicates that the time At'^ and the time t'a^^ are 
affected by the position of the boundary due to the boundary-dependent term a'^^^^-^ in 



Eq. (20) and (21). In the first order of perturbation w.r.t. Gc-, all deviations from the 
SFA saddle point solutions are proportional to Gc'- AO'^ ~ 0{Gc), Acf)'^ ~ 0{Gc) and 
At'^ ~ 0{Gc)- In our analysis we shall consistently keep only terms of the order ~ 0{Gc) 
and therefore, only the SFA saddle points can enter the argument of Gc{p, T; a, t') and the 



SFA action. Indeed using Eq. (21) we obtain: 

Gc(p, T; a, t') = Gc (p, T; rf ), tf )) + O(G^) (23) 
S'^\p, T; a, t') = 5^^^ (p, T; rf ), tf )) + O(G^) (24) 

However, the corrected saddle point solution for time will contribute to the pre-exponential 

[e'4 

dt' 



factor P{6' ,t'), since ^^%lf '* 7^ 0. It is straightforward to show that = 6^{t') and 



t>v{t'), where 6^{t') and 0t,(t') describe the direction of electron velocity at time t'. 



C. Integration over the surface sphere 

Because of the large value of the action, only a vicinity of saddle points contributes to the 
integral. We make Taylor expansion of Gc around points 6'j^\ (f)s^\ t'a^^ (only the saddle 
point in time is affected by the boundary term avp{t')) up to quadratic terms: 

Gcip, T; a, t') = Gc (p, T; rf ), tf )) + (a - rf )) • VGc (p, rf ), t^^")) 

+ {f - tf )) d,Gc (p, T; rf ), )) + i (f - e^fd^Gc (p, T; rf ), tf )) ^^^^ 

The mixed derivative VdfGc [p-, T; r'i°\ t'i^'' j oc C(G^) is omitted from Eq. (25), as higher 
order correction, since it is multiplied to t'—t'a^^ oc 0{Gc)- The term involving second deriva- 
tives w.r.t. spatial coordinates on the surface sphere ^ ^a — rs^"* j AGc (^,T;r'j^\ta'^^^ 
is equal to zero for Coulomb potential, since AU{r) = 6{r) and the argument of U in 
Gc (^P,T;r's^\t'a^^^ is trajectory starting at the surface and propagating outside of the 
sphere; this trajectory never reaches the origin. Note that VGc = ^P? where Ap is the 
shift of canonical momentum due to electron interaction with the long range potential of 
the core: 

Ap(t', T) ^VGc = - r dr .^.r^' [r' + F d( Vp(C)l (26) 

A' \\^ + Jt' ^Cvp(C)|| L Jt' 
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where U' represents derivative of U w.r.t. its argument. It is convenient to rewrite the time 
derivative of Gcijp^ T; a, t') as P 



d,Gc (p, T; rf ), tf )) = -Vp (tf )) ■ VGc (p, T; rf ), )) + f/(a) ^^^^ 

= -Vp-Ap (tf),T)+f/(a) 



Substituting Eq. (25) into the expression for ap(T), and evaluating the integral over 0' and 



6' exactly (see Appendix B 1 ) we obtain: 

ap(T) = Ar,„(-.)^(-l)-^,K«)^e^''"("'^^"'° T rft' e't /"'^^'^p«+-^(*'-*o)/2 

V27r (28) 

\fpc(t')/ 

where (py{t') is the tunnelling angle: 

and p'^ = p — Ap(i(:', T) is the Coulomb-shifted momentum at time t' corresponding to 
asymptotic momentum p registered at the detector at time T: 

Ap(t',r) = VGc(p,T;rf),t') (30) 

Note that the tunnelling angle is complex. It signifies sensitivity of strong-field ionisation 
to the sense of rotation of the electron in the initial state [HI El] . 



D. Integration over time 



We are now left with the integral over t'. We will use the saddle point method to evaluate 
this integral. The saddle point equation for t' is: 



dS 9^SFA QQ 



df 



dt' 



+ 



dt 



/(0)^ 



c 



+ av'(t')=Q 



(31) 



where the last term, as discussed in subsection HIE, comes from jg {avp{t')). To solve this 



equation, we expand the derivative of the SPA action in Eq. (31) up to quadratic terms 



w.r.t. At'i^^ and take into account that dt'S^^-^ (ta^^) = 0, yielding: 



d,Gc (tf ) 



d,Gc (tf d,Gc (tf ) 



df,Sr^{t'P]+df,Gc{t'} 



(32) 
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Here we have used that t'i"^ — t'i^^ = —ia/n <^ t'}'^' [TO]. We have omitted terms of order of 



/(O) 



Gc in the denominator in the last term in Eq. (32), since the terms of the first order of Gc 



are already included in the denominator. Taking into account Eq. (27) and 



(33) 



we obtain for the time: 



-Vp (tf^)-Ap [e\T)+Uia) 



E t 



Vp I t 



/(O) 



(34) 



Eq. (23) and (24) suggest that Coulomb corrections to the ionisation time do not affect the 



exponent of the ionisation amplitude, but they contribute to the prefactor, further modifying 
the tunnelling angle: 



tan0S (tr + Atf)) 



p, - Ap, {t'P + Ae\ r) + Ay {tT + Atf) 

p. - Ap^ (tf ) + Atf t) + (tf ) + Atf ) 



(35) 



Up to first order terms w.r.t. to Gc, cip{T) is: 



ap(T) = iKa^(f^i{a)Nim 



-t /T(0) drvliT)+iKh'}''^/2+iGc(p,T;r'}''\tT'> 



-<'"''hr[^;^]je {av, (if))) 



Vr,c to: 



(36) 



E. Boundary matching 

We now consider the elimination of boundary-dependence in the results for transition 
amplitude. In the long pulse, due to cylindrical symmetry of the problem, the result does 
not depend on position of the detector in the polarisation plane x, y. Thus, without the 
loss of generality we will consider electron registered at the detector placed in the positive 
direction of x-axis, i.e. the electron momentum at the detector py = 0. 
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1. Complex momentum Ap [ta^\T\ at the boundary 



To perform boundary matching in Eq. (34) and Eq. (35), we need to evaluate the mo- 
mentum Ap (ta'^\T^ at the boundary: 



Ap(t'i°),T) 



/(O) 



/(O) 



(37) 



Note that Ap (tT\T^ is a function of final momentum p. In this section we consider only 
P = Popt; corresponding to the momentum at which the probability is maximal, since it is 
sufficient to calculate the ionisation rates. The photoelectron spectra will be considered in 
our subsequent publication [T^. In the polarisation plane the optimal momentum popt = 
(Popt COS 0p, Popt sin0p) is given by radial momentum 




Popt 



for any angle (pp. The parameter < Co ^ 1 satisfies the equation 



tanh 



(38) 



eg +7^ 



pin [21]. Note that Co ^ for 7 < 1, and Co ^ 1- l/ln7 for 7 > 1 [14j. An alternative 
expression for popt is 



Popt = ^0 



sinh CUT' 

w 



/(O) 



(39) 



UJT, 



where rf) = c>[tf)],is the imaginary part of the saddle point solution for time, also known 
as the "tunnelling time". The advantage of the second expression is that it provides a 
compact connection between the optimal momentum and the tunnelling time, however, one 
has to keep in mind that in circular field t^^^'' depends on final radial momentum p [HI [TU El] : 



OJT, 



cosh""^ T], ?7(p) 



p_ 



+ 7+1 



(40) 



and thus, in Eq. (39), t^^^^ depends on popt itself. 

Since the time t'i^^ is complex, the momentum Ap ^t'i^^T^ will also be complex: 



zAp''{a) 



Apy {tf\T) = Ap'^{a 
Ap. {tT\T) = Ap'M) + ^^p:{a 



(41) 
(42) 
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After some algebra (see Appendix A 1 ) we obtain: 



Ap''{a) ^O( — ]^0 
^ \Ka J 

Ap'^ia) =Ap'^, Ap'yia) = Ap'^ 



(43) 
(44) 
(45) 



where the boundary independent momentum is 



Ap' = - 



and the coordinate re^"*, known as the coordinate of exit from the tunnelhng barrier, is 
defined as 



drVU r'i°) + 



.'{(>) 



rfCvp.,,(C) 



(46) 



.'(0) 



1(0) 



'(0)l 



C^CVp„p,(C) 



Finally (Appendix A 1 ) 



(47) 



(48) 



2. Boundary matching for At'a'^^ in Eq. (34) and for the tunnelling angle tan cp'^ (t'a^^ + At'} 



AO) 



AO) 



in Eq. (58) 



Substituting Eqs. (pi), M^, (|45|), (|46l), (pi) for Ap ( tf \ T ) from the previous section 



into Eq. (34) and taking into account that in our geometry 



E,(tf)) =4, 4 = Eocosh 



LOT- 



'(0) 



(tf )) = iE'^, E'^ = -Eq sinh wr. 



/(O) 



—Aq sinhcur^ 



'(0) 



Popt - Aq cosh UT^ 



/(O) _ «0 



'(0) 



sinh UT, 



/(O) 



UT cosh UT, 



/(O) 



where oq = Eq/cu is the electron oscillation amplitude, yielding 



-'(0) 



E {tT^) ■ Vp„p, [t'P] = zpopt^osinhwr.^ ' = iVyPo^ti^ 



we obtain: 



53[At'i 



'(0)1 



PoptUJ 



(49) 
(50) 
(51) 

(52) 



(53) 



(54) 
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Since U{a)/Vy ~ Apy ~ O (^) — )■ (see Appendix Al), Ap" = we obtain boundary- 
independent correction to real ionisation time 3fJ[At'i°^] = 3?[Ata°''] and imaginary ionisation 
time 3[Atf = '^[AtT]: 



Ap'^v'^ Ap'^ sinh ut^^^^ — ujt cosh ojtI^'^^ 

V''PaptUJ Eq 



sinh^ wr,'*-*^^ 



3?[At: 



'(0)1 



Ap; _ Ap'y ur[ 



(0) 



-^0 sinh wr,'*^'^'' 



(55) 
(56) 



where the subscript "s" denotes that the resuhs for corrections to SEA saddle point t'lf'^ are 
now independent of the boundary r' = a. Thus we can write the saddle point as: 



(57) 



Matching for the tunnelling angle is now trivial, since all variables entering Eq. (35) are now 
proved to be boundary- independent: 



tan<^^ (tr + Atr) 



Vy {t?^) - Apy - At'PEy 

V. {t'P) - Ap^ - AtT^E^ 



(58) 



3. Boundary matching of remaining terms in Eq. ( 36 ) 



We first establish the connection (in Appendix A 2): 



(59) 



Next, we consider matching of the EVA phase to the bound wave-function. We follow the 
approach used in [TU] for linearly polarised field. 



B{a) = na (p,,e{a)e V Jj^, i^avp [t^ >)) 



(60) 



The asymptotic bounded wavefunction in a Coulomb potential is: 

Eurthermore, at r' = a, we can write (na)'^^'^ as 

[Kap' = e •'1/'= X = e ^ " ^ = e '« ^ ^ 



(61) 



(62) 



15 



where t^^^ = t's^^ — i/k"^ and ta""* = t'j'^'' — ia/n. The second equahty follows from the fact 
that between t'i"'' and ta^\ the velocity of the electron remains almost constant, while the 
third equality holds because finally we will be using the modulus of the vector and under 
the approximation: 



'(0) 



'{0) 



(63) 



The term given by Eq. (62 ) can now be matched with the Coulomb phase term Gc yP, T; r's^\ta^^ 

ft?' / r \ r'-" I r \ 

Gc(rTKe\T)=J^ dru{r'r + l^JC.M)=l drf/ dCv,(C)) (64) 

to yield after using the large argument approximation for ji (^aVp (t's^^^^: 

which ensures boundary matching for all transition rates/amplitudes. 

After boundary matching the final expression for the ionisation amplitude is independent 



of a: 



-5 Ao) drvl{r)+i^tT^ -i/J(o) dr C/(/;(o) dCvp(C)) 

— — g g '■hi V *S / 



S" t'P 



(66) 



-I 



Pz 



'(1) 



Here t's^^ is given by Eq. (57). Eq. (66) corresponds to "single" ionisation event (ionisation 



amplitude formed after one laser cycle), since only one saddle point is included. 



F. Ionisation rate 

To calculate the ionisation rate we integrate the ionisation amplitude (corresponding to 
single ionisation event) over all momenta using the saddle point method and divide by period 
of the laser field: 

(67) 



27r 



w = ^ I dp \ap{T)f 
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In Appendix B 4 we evaluate this integral exactly, but consider here the rate around the 



PpC^^p and the integral Eq. (67) can be calculated using 



optimal momentum. As it follows from Eq. (66), the ionisation amplitude ap{T) can be 
written in the form: ap{T) 
saddle point method. 

The saddle point equation 



Vp253[Fp] = 2Vp3[5f ^] + 2Vp [53[Fp] - ^]] = 



(68) 



can again be solved iteratively, since the second term is small by construction. The optimal 
momentum in SEA solves the equation: 



OrcSFAj 







(69) 



and is given by Eqs. (38), (39) Since the correction to popt are obtained from the Eq. (68), 
they will contribute to the ionisation rate in the second order w.r.t. Gc- We keep only terms 
of the first order w.r.t. Gc and therefore these corrections are irrelevant and the saddle 
point for the momentum integral in the ionisation rate is given by the optimal momentum 



Eq. (38). We neglect here small corrections arising from substituting the pre-exponential 



factor S'^ (tf ^) by S'^ (tf ^) in Eq. (|66|). Finally, using saddle point method for radial 
integral and taking into account that integration over (pp yields 27r we obtain the expression 
for the ionisation rate: 



opt 



7 



TT 



y;^V9[^"(Pp,opt)] 



1^ ^2Wci+2Wc2 




-2mS 



(70) 



where 



5>[^"(Popt)] 



2Co^ + 7^(l + Co^; 



c^(l-Co)V(a + 72)(l + 7^ 



Wci 

Wc2 



U 



/(O) 



^Cvp(C) 



'{0) 



U 



'(0) 



^Cvp(C) 



(71) 

(72) 
(73) 



Wci is well-known adiabatic Coulomb correction, evaluated under the barrier along the 



optimal trajectory [Ml [15]. Analysis of Eq. (70) shows that non-adiabatic Coulomb effects 
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modify the ionisation dynamics in several ways. New effects arising from our anafysis include: 
modification of (i) ionisation times, (ii) initial conditions for electron continuum dynamics, 
(iii) the "tunnelling angle" . We discuss these Coulomb effects in detail in the next section. 
We show that Coulomb effects modify (i) calibration of the attoclock in the angular streaking 
method, and (ii) the ratio of ionisation rates from p~ and orbitals obtained for short range 
potentials in [H] . The photoelectron spectra will be considered in our subsequent publication 
[19] , where we will include the effects of Wc2, the result of interaction of long-range potential 
with the electron in the continuum. 

IV. PHYSICAL PICTURE OF IONISATION IN LONG RANGE POTENTIALS 

In circularly polarised fields, the electron liberated at different times will be "directed" 
by the laser field into different angles. This idea is called 'angular streaking' and the corre- 
sponding "time-to-angle" mapping is unique for nearly single-cycle pulses with stable carrier- 
envelope phase, underlying the idea of the atto-clock P-H] . The angular streaking principle 
makes single and double ionisation in circularly polarised strong laser fields a sensitive probe 
of the attosecond dynamics P-[6]. 

However, reconstruction of this dynamics requires the calibration of the attoclock, i.e. 
establishing the mapping between the direction of the laser polarisation vector at the time 
of ionisation and the direction of the electron momentum at the detector. When one strives 
to achieve the accuracy of, say, 10 attoseconds, using 800 nm carrier as a clock, one needs 
to know this mapping with accuracy of about one degree. 

Simple analytical calibration can be made if one neglects the electron interaction with 
the the long-range core potential during and after ionisation. For short range potentials the 
mapping is illustrated in Fig. [TJ For the laser field defined as: 

E{t) = Eo{- sm{ut) X + cos(wt) y) (74) 

the connection between the real part of the ionisation time and the observation angle is 
[12 [21]: 

wtf^ = u^tf^] = 0p + 2n{r - 1), r G N (75) 

The detector placed at positive direction of x-axis will detect the electron liberated at 
= 0, i.e. when the laser field E(t) = Eoy is pointing towards the positive direction of 
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y-axis. The electron exits the barrier in the negative direction of y-axis, corresponding to 
the angle —tx/2. The velocity at the exit Vy j = 0, {^i^^^ = Popt — ^o? and [t'i^^^ 
tends to zero in tunnelling limit (7 ^ 1): Vx [t'i^''^ = a/2/p7/6. Thus, the angle between 
the direction of the field at the moment of ionisation and the electron momentum at the 
detector is 7r/2. 

How this mapping is affected when the interaction with the long-range core potential is 
taken into account? 



A. Coulomb correction to the ionisation time, initial electron velocity and the 
calibration of the attoclock 

Even in the tunnelling limit, our analysis shows that due to the effects of the long-range 
potential the electron has non-zero velocity {—Apy) in the negative direction of the y-axis 
when the field is pointing in the positive y-direction, i.e. at t = in our notations. 

This is by no means surprising and the corresponding velocity has a very simple explana- 
tion: it is required to overcome the attraction of the Coulomb potential, which the electron 
will experience all the way towards the detector. Had the electron been born with zero 
velocity in long range potential, it would have never reached the detector placed in positive 
direction of a;-axis. One expects the same result within the adiabatic tunnelling picture. The 
question is: is the magnitude of Apy consistent with the adiabatic ionisation model, which 
would suggest that the electron was liberated slightly before t = but with zero velocity? 

To answer this question, we need to analyse the changes in the ionisation time due to 
the effects of the long-range potential. The corrections to ionisation times associated with 



electron interaction with the long-range potential are given by Eqs. (55) and (56). The shift 
of the saddle point in time 3ft[At'i°^] corresponds to the shift in the direction of the force of 
the electric field -E(t) from -7r/2 to -7r/2 + uj^[/\tT\ 

Let us first discuss the initial conditions for the electron continuum dynamics in the 
tunnelling limit 7^1. In this limit, the momentum shift is accumulated along the electron 
trajectory: 



z.(t-tf)) + |^(t-tf))^ 



(76) 

The tunnelling trajectory along y-axis can be further rewritten as: ']f'^^{t) = —Ip/E^ — 
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^Ef) (t — t'/°^j . Hence, we are only left with the real part where y*""^ (^tf^^^ = ^Ip/Eo is the 
coordinate of the exit point in the tunnelling limit. Taking into account that U = —Q/{—y), 



VU = —Qy/y"^ and substituting this trajectory into the expression for Ap, Eq. (46), we 
obtain: 

mPy] = -Q fl = -'-^QEo, (77) 

'(0)i 



/j/(0) (ytunp j3/2 

It is easy to see that Eq. (56) in the tunnelling limit yields 9f?[At'i"''] = —Ap' /Eo, thus we 



obtain from Eq. (77): 



.3/2 - > C^S) 



From Eq. (78), we find that the correction to the ionisation time 9fJ[At'i°'*] is negative, the 
electron is born before E(t) points down, and the Coulomb corrected angle —7r/2+ui^[At'j^''] 
has a larger negative value. At this (earlier) ionisation time the electron velocity is smaller 
than at tf^^ and in the tunnelling limit: 

V, = Popt - Ao cos ) + u;3ft[Atf )]) - Ap', ^ popt - Ao + O(G^) ^ 7 v^/6 + O(G^) 

(79) 

-Apy - Ao sin (^wtf ^ + w3ft[Atf )]) = -Apy - Aocu^At'^] ^ + 0(G'^) (80) 



Vy 



Thus, in the tunnelling limit 7 — )■ 0, the electron velocity indeed tends to zero at the exit 
from the barrier. The effect of Coulomb potential is reduced to the modification of the angle 
between the direction of the laser field at the moment of exit E j and the direction of 
final electron momentum p, registered at the detector. For short range potentials this angle 
is 7r/2, for long range potential this angle is larger; in tunnelling limit it is n/2 + ujlp , see 
Fig.[Tj 

However, most of the experiments are currently performed in the regime of non-adiabatic 
ionisation, when the Keldysh parameter 7 is not that small. In this regime the exit velocities: 

V. = Popt - Ao cos (a;3fi[Atf )]) - Ap, (81) 
Vy = -Apy - Ao sin {uj^[Atf^]) (82) 

become significant already for small 7. The longitudinal electron velocity v\\ along the 
direction of the field and the electron transverse velocity v± orthogonal to the field are also 
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7r/2 - Act 



Vy^in sinAa 



I, in 



(a) Short Range Potential 



Vx,in sin Aa + Vy^in cos Aa 
(b)Long Range Potential 



FIG. 1. Kinematics of electron tunnelling through rotating barrier. Right circularly polarised 
laser field E creates the tunnelling barrier rotating counter-clockwise, (a) Short range potential: 
The electron observed at the detector placed along the j;-axis, exits the barrier along the negative 
direction of y-axis at ang le a(0) = -7r/2. (b) Long range potential: The electron observed at the 
detector placed along the x-axis, exits the barrier at the angle a^^'^ = —tt/2 — Aa, Aa = uoAt 
and a6^^ < 0. 



'(0) 



non-zero (Fig. [2]). The longitudinal and transverse velocities are obtained from Eqs. (81) 



and (82) (Aa 



/(O) 



v± = Vx COS (Aa) — Vy sin (Aa) 
v\i = Vx sin (Aa) + Vy cos (Aa) 



(83) 
(84) 



Ignoring non-zero initial velocity of the electron will lead to errors in the two-step recon- 
struction of time delays in angular streaking method. 



B. Coulomb correction to the electron "tunnelling angle" 

The complex tunnelling angle characterises the direction of the electron velocity at the 
complex ionisation time t's^^: tan0t,(t'i^^) = • The ionisation rate is proportional 

to the imaginary part of the tunnelling angle w oc e2m,3[(/>^(ts )]^ where m is the magnetic 
quantum number. In case of spherically symmetric initial state (s-state) m = and the 
ionisation rate does not depend on the tunnelling angle, because the electron density in 
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Field frequency O), eV 



FIG. 2. Initial velocity corresponding to the center of the velocity distribution: v± (green dashed) 
(Eq. ISSl) and vu (blue dot-dashed) (Eq. Isil) vs frequency for Eq = 3 x 10^° V/m (^o = 0.06 a.u. 



and I = 2.6 x lO"*^^ W/cm^) and Ip = 14 eV. (red) shows the result for short-range potentials 
(SFA velocity). 

the initial state is the same in all directions. For p-states however, the direction of electron 
tunnelling, defined by the tunnelling angle, becomes important. In particular, it leads to the 
sensitivity of ionisation to the sense of rotation of the electron in the initial state. For short 
range potentials this effect was predicted and analysed in In this section we discuss the 
non-adiabatic Coulomb corrections to the tunnelling angle and show how the results of [8] 
are affected by the electron interaction with the long-range core potential. 
The tunnelling angle in case of short-range potentials is: 



Py - Aosin (utT^ 
tan 0„ (tf )) = y—L (85) 

Px — Aq cos f Uts 



The Coulomb potential leads to the two equally important effects: (i) the modification of 
the complex ionisation time {t's'^^ + At'i^"* in the long-range potential vs. just t'}^^ in the 
short range potential), and (ii) the momentum shift due to the deceleration of the electron 



by the long-range potential of the core (see derivation in subsection HIE): 



vje^)-Apy~Ae^Ey 



(tf ) j - Apx - Atf^Ex 
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In this section we focus on the imaginary part of the complex tunnelhng angle (pHt'g) = 
tan"^ {x + iy), since it contributes to the ionisation probability. The imaginary part of (f>y{t'^) 
can be cast in the form: 



^^M)] = -\ In [{l-x'-y'Y + Ax') + ^ In ((1 + yf + x') 



(87) 



Note that the real part x ^ 0{Gc) is of the first order with respect to long-range potential 
and therefore the terms have to be omitted. The ratio between ionisation rates for p~ 
and p"^ orbit als is 



,4S[</,S(t'i'')] 



l + y_ 



Finally, 



y 



< + APx</bopt tanh ut] 



(88) 
(89) 

(90) 



Fig. |3] shows how the non-adiabatic Coulomb effects change the ratio between the ionisation 
rates for the p'^ and p' orbitals. Modifications come solely from the alteration of the 
tunnelling angle. The non-adiabatic Coulomb corrections {Wci and Wc2) do not contribute 
to the ratio of the ionisation rates, as has been also noted in [8]. The decrease in the p~ /p'^ 
ratio at high frequencies in long-range potentials is consistent with the opposite propensity 
rules in one photon ionisation, where p~^ is preferred over p~ for right circularly polarised 
fields. 



V. CONCLUSION 

We have evaluated strong field ionisation rates and amplitudes for circular fields taking 
into account non-adiabatic barrier dynamics of a Coulomb potential using the recently de- 
veloped ARM technique. The ionisation rates for atoms in arbitrary potentials in circular 
fields for long range potentials have been derived rigorously, extending the work of [HI (TD] 
and [m |22]. Along with the Eikonal-Volkov approximation to the free electron state, the 
ARM approach allows for accurate and rigorous analysis of ionisation in strong fields, con- 
sistently including the Coulomb effects both during and after ionisation. It should be noted 
that in the current implementation of the ARM we have included Coulomb effects in first 
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^ 6 



■BS-2011 (short-range) 
■ARM (long-range) 




2 4 6 

Field frequency (O, eV 



FIG. 3. The ratio of ionisation rates from p~ and p"*" orbitals for Ne atom {Ip = 21.5645 eV) and 
Eo = 7.7 X 10-'^*' V/m {Eq = 0.15 a.u. and / = 1.6 x 10^^ W/cm^): Wp-/Wp+ for right circularly 
polarised field: short range potential (red) and long range potential (blue dashed). 

order perturbation to the action. It limits the applicability of current implementation to 
the region of moderate 7. The simplest "post-mortem" validity check can be performed by 



computing Ap' (Eq. 46 ) and comparing it to the SPA velocities. The momentum shifts Ap' 



should not exceed the SEA velocities. 
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Appendix A: Supplementary information for boundary matching 



1. Complex momentum shifts at the boundary 

The goal of this section is to calculate momentum shift at the matching point a: 
Ap(a) = -^^^rfrVf/(rf) + ^^^rfCvp(C)) 



(Al) 
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and show that it does not depend on the position of the boundary under the matching 
conditions. We first spht the integral into two parts: 

r.5R[t'i°'] 



Ap(a) = - 



,'(0) 



drVU{r'P+ / , rfCvp(C) 



.'(0) 



'(0), 



1(0) 



Physically, these two parts can be interpreted as accumulated before 



and after: 



ApW(a 



Ap(2)(a 



'(0) 



1(0) 



'(0), 



dTVu(r'P+ / rfCvp(C) 



'(0) 



(A2) 



(A3) 



(A4) 



the tunnel exit defined as the coordinate at the time ^[t's 



(0)i 



.'(0) 



1(0) 



'(0), 



(A5) 



The second part, Ap*^^)(a), does not depend on the boundary. In the following we will show 
that the first part Ap*^^)(a) is negligible under the matching condition na ^ 1. We first note 
that Apl^\a) is purely imaginary, while Ap^x\a) is purely real. In the same geometry that 
we have used in the main text of the paper t'a^^ = iTa^\ and the complex under-the-barrier 
trajectory is R = r + ip: 



-ao 



p = ao 



cosh (j)'}^^ — cosh < 



sinh d)'}^^ — sinh < 



(A6) 
(AT) 



where (pf^'^ = ujt[^^'' and = cj^ and ^ is imaginary integration time variable. The Coulomb 
potential takes the form (details of the analytical continuation of the Coulomb potential to 
the complex plane will be addressed in our subsequent publication [IS]): 

1 



UK) 



(AJ 



Purely imaginary Ap"y\a) is 



Ap«(a 



(A9) 
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Purely real Api\a) is 



fO 

.(1) 



= TT^^^- (AlO) 



1(0) 



J.2 _ p2^3/2 



and in both cases, (pa = ^t'o^^ ■ Also, since for the optimal trajectory r ^ p, 

Ap«(a) ^ r '-f. (All) 

As /O = at the tunnel entrance ^0 = 0s = wts''^ j and p = at the tunnel exit (0 = 0), the 
integral is accumulated in the vicinity of Ta^\ We make linear expansion of the integrand 
around this point: 

where C is a numerical factor: 

C = / 7-^- (A13) 



.0 (^ + 1)^' 

So far we have considered Ap(a) defined through its outer region value. We can also 
estimate Ap(a) using its inner region value. The inner region value of Ap(a) can be calcu- 
lated using static approximation (or short time propagation), since the time interval from 
t's^^ to t'i^'' is very small. It is convenient to estimate Apy{a) by evaluating its inner region 
value. In static field, the momentum in the inner region p"(a) is defined through the energy 
conservation: 



.p- E^a (A14) 

la 

Thus, j}^{a) = -i^2{Ip - Eoa - 1/a) ~ -i^2{Ip - Eoa)il + l/(2a(/p - ^o^))), yielding 
p™(a) = ~iK{a) — i/K,{a)a. The first term is the SFA velocity at the boundary K,{a) = 



a/ 2(/p — Eoa) ^ K, the second term is the respective correction associated with Coulomb 
effects. Thus, Ap™(a) ^ 0{1/K,a). The vanishingly small value of the correction at the 
boundary is not surprising, since the boundary is placed in the region, where the Coulomb 
modification to the barrier is already very small (see Fig. |4]). 

2. Additional expressions for boundary matching 

We derive here the relation: 

Ji (aV e-'"'-^^ = Je {av, (tf ))) (A15) 
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FIG. 4. Cartoon picture of the matching philosophy. As indicated by the two red arrows the 
imaginary momentum Ap^^\a) is much smaller in the matching region than in the inner region. 

Since the saddle point t'i^^ is close to the SFA saddle point t'j^\ we know that the argument 
of je is of the order of ^ 1. So using the large argument approximation for spherical 
Bessel function, and expanding Vpc(t') up to first order in Ap, we get: 

W'O) = Je K (A16) 

It can be shown that: 

Since the inner region should be treated in quasistatic approximation, the second term is 
vanishingly small. 

Analogously to boundary matching approximation made in [20] we obtain: 

{tT) - m (A18) 

Taking into account that by definition 



and 




(A19) 



(A20) 
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au 



AtfV(C^ + 7^)(l + 7^) 



(A17) 



we obtain Eq. (A15). It must be noted that it was because of Eq. (A15) that je{cLVp{t')) 
was used in Eq. (31) and not ji {aVpc{t')). 



Appendix B: Prequency-domain approach 



1. The wave-function 



In our frequency-domain approach we start the analysis with the expression for the wave- 
function in the coordinate representation ipir^t) = {r\ipont(t)) , where \4'out(t)) is given by 
Eq. (g: 

i^out{r,t)=zl'dt' I dr' I rfr"(r|f/B(t,t')|r')(r1^-(a)|r")(r"|^in(t')) (Bl) 

Taking into account the explicit form of the Bloch operator in coordinate representation 
Eq. ^ we can rewrite Eq. (Bl) as follows: 



dt' j dr'Gsir, t; r' , t')6{r' - a)B{a, 9', 0', t') (B2) 



After following the arguments in subsection III A, we can approximate the boundary term 



B{a,6',(j)',t') as in Eq. (18). This follows from the foresight that when we use the saddle 
point for the time integral to isolate tunnelling dynamics, we end up with studying the 
dynamics of the wavefunction around the pole v{t'g) = iK in the momentum space. This 
corresponds to prominent contribution only from the asymptotic part of the wavefunction 
in the region Kr ^ 1 in co-ordinate space. 



Using Eq. (18) and Eq. (25) and evaluating the Delta function over r', we now have for 



the wavefunction ilj{r,t): 

^(r,t) = dt' dk dO' / t/0'e*(^>^W-^>^(*')-^)e-5i;'^-"^M^,,(a) 



where Nim = y [e+\m\y. ^ ~ a(sin 6' cos 0' x + sin 6' sin 0' y + cos 6' z) . We also take 
A(t) = — 74o(cosa;tx + sinwty)- 

The point r'i^-* is defined in spherical co-ordinates as (^a,9'j^\(j)'s^^^, where = 9^(1:') 
and 0s = (f)v{t'), which gives: 

rf) = a^^^ f dr^r^ir) (B4) 
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The approximation follows from the fact that the saddle point for t' will be quite close to 
the SFA saddle point ts°^ and hence by defining t' — t'j'^^ = — rnr, we can redefine saddle 
point r'i°^ as a classical trajectory, t' here corresponds to zeroth order correction in the SFA 
saddle point when the electron is propagated from a finite boundary instead of the origin 
(thus the saddle point for S^^^ + avp{t')). 



Following subsection IIIC , the resulting surface integral is: 

/•TT /»27r 

I^,= de' sine' rf0'e-^^(*')-^Pf (cos^')e'""^'e*^''"' (B5) 
Jo Jo 

The term e*^'' ** comes from the Taylor expansion of the Coulomb phase Gc about the saddle 
point co-ordinate (^a,9'}^\(f)'s^''^. Since the gradient of Gc is identified as the momentum 
shift, we can see this term as the contribution of the long-range potential to propagation from 
the finite boundary r' = a. Including this shift by rewriting the shifted kinetic momentum 
as Vkc(it) = Vk + A(t) — Ak, the integral over 0' is evaluated to: 



^0 



= 27re^™*S(*') (avkc(t') sin^') 6"^'^'=^=°^^' 

The superscript "c" denotes that we are calculating surface integral over Coulomb shifted 
momentum and J„(-2) is the nth-order Bessel function of the first kind. 
The f2'-integral now is: 

/•TT 

I^, = 27r(-i)'"e'™'^sa') / de' Jm {avi,c{t') sinO') {cos 9') e-^-K-^^^' gin^' (B6) 

Jo '' 

We depart here from the method used in [TU] of approximating the 6''-integral around a given 
angle according to the direction of polarisation (there, 0' ~ tt was a reasonable approxima- 
tion, and here Q' ~ 7r/2). But with the & ~ 7r/2 approximation, not only do we lose accuracy 
in our result, but the small-argument approximation would not be valid for Jm(a6 sin 6''). 
But we have used Q' ^ -n j2 for the Coulomb correction, as deviation from a planar trajectory 
here is suppressed exponentially J23j. Hence we perform an exact analysis, noting that the 
above integral has an analytic expression from by using a similar integral on product of 
Bessel functions and Gegenbauer polynomial from [25], which finally gives us: 

l^, = 47r(-^)^(-l)-e^'"*S(*')p- (-^] {av^^it')) (B7) 
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Substituting this result into Eq. (B7) and using Appendix A 2, we get the wavefunction as: 



A dr U{rL (r ;r,k,t))+iGc (k,T;r'i°' ,t'j p 



(B8) 



w ji (afk(t')) 



Eq. (B8) is an exact expression from ARM under the PPT approximation. 



2. lonisation Rate 

In order to calculate the ionisation rate, we need to know the radial current density, 
jp(r, t), defined as 



(B9) 



Following the procedure of fl^ , but noting the changes due to the presence of the Coulomb 
phase term, we can get the familiar expression: 



oo „ 

w{S,uj) = 2ttJ2 dk\Fn{k'',uj)\^6 

n>no 



r 



nu 



(BIO) 



with 



F„(k^a;) 



2^ 



27r 



dt' F{k'',t')e 



(Bll) 



The Coulomb phase term is the main difference from the result for short range potential. 
3. Derivation of F„(k, 



Unlike the result for short-range potentials [H], we now have an additional term in 
the exponential oscillations due to je{av]i{t')), along with the Coulomb corrections. Apart 
from the modified. Coulomb-shifted momentum that is a new result from this analysis, the 
Coulomb term in the action also includes motion after ionisation, introducing a modification 



of the result in [8], [HI [221 [26] . As discussed in subsection III E , we know that the saddle 
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point in time would be such that v {i!a^^^ ~ ±m, and same as done there, can use the 
asymptotic condition for large argument {na ^ 1) on the spherical Bessel function: 

The two terms correspond to contributions from the diametrically opposite points on the 
boundary surface a, from where we propagate the electron outwards. The point farther 
from the detector by a distance of 2a compared to the point nearer causes an additional 
exponential decay for propagation from the former. Such a term didn't appear in [TU], as 
there saddle-point analysis on k-integral was used, thus isolating the electron field to one 
particular trajectory, corresponding to a classical particle motion rather than field evolution. 
Not using saddle-point in our case would naturally lead to inteference effects between the 
contribution from the two points, but under the given condition (/ta ^ 1) those effects will 
be exponentially small. This way, an interference will be produced on every point throughout 
every circular disc for different 9 on the sphere r' = a. The contribution of each is weighed 
by the momentum distribution, encoded in e«"^</'S(*')p™ xhe maximum contribution 

comes from the region around the saddle-point, which effectively considers the electron as 
a particle. However, since our analysis is exact, the contribution from momenta about the 
classical are also included in the above result, as well as taking into account the case for 
non-zero perpendicular momentum {k^ 7^ 0). 

The saddle point corresponding to the boundary-dependent action S^^^ = S^^^ + aVp{t') 
can be derived after Taylor expansion about the SFA saddle point t'j'^^: 



tf)=tf)-^^ (B13) 
After modifying the SFA saddle point t's^^ through the change in t'a'^^ due to the Coulomb 



phase term, as discussed in section III , we get the final expression for the n-photon transition 
amplitude, to first order in a: 



^'»(k.-) = |^(->)'(-ir 



AT -i-^uyi^s I pm I '"z 1 im<j>%[ts 



(B14) 
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After boundary-matching (subsection 



HIE): 



27r 



^7 ^-»5o(tf))+fm0g(f'i^>) 



(B15) 



Since we are interested in |F„(k, only, we get: 



UJ-i2i+l{i- \m\)\ (1 - (-1)^^"")' 
167r3 (£+|m|)! 



V]^c t 



.'(1) 



(B16) 



-2mS [<^S (t^(^')] -2^ (r, cosh"! - V^) 21^01 



+2Wc2 



For short-range potentials {U = 0) the above result matches equation (17) in [8J precisely. 

We see another advantage of the ARM here: we now don't have a complicated radial 
r'-integral and the corresponding higher-order pole in momentum space representation of 
the wavefunction. The upshot of analysis in short-range potentials |H] was that the pole in 
the momentum space representation of the wavefunction was cancelled with the zero in the 
momentum integral at the same point v{t'g) = in. However, for wavefunctions corresponding 
to long range potentials, we would have had a Q/n + 1 order pole in the momentum space, 
leaving a. Q/k order pole in the final momentum integral. Using ARM, the Bloch operator 
isolates the wavefunction at the boundary r' = a through a delta function, making that 
integral straightforward, thus bypassing the pole encountered if the integral was performed 
over the whole radial domain. At the same time we also get a more robust result, taking 
into account Coulomb correction for the ionisation rate both during and after ionisation. 
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4. N-Photon lonisation Rate 



The n-photon ionisation rate is 
Wni£,uj) = 2n J (ik |F„(k, cj)!^ 5 



87r2 {i+\m\)\ 



^2Wci+2Wc2§ 




2ir poo 

d(f)k / 

^0 



(B17) 



nu 



2 \ \ T 

Using the Delta function, the integral over kp is easily done by substituting kp = \/k\ — k\^ 
where k\ = 2nu - k'^ (^1 + We modify the definition of C = - l), used in [22] to 
include contribution from trajectory perpendicular to plane of polarisation to give: 



eff 



2nf 
n 



1 



(B18) 



where 2n'^uj 



^eff 



5 "-eff 



+ k1 and 7eff = Hesl^o as defined before. The 



corresponding values for different functions of k appearing above are as follows: 



1 + 7e^ff 



1 - 



'1- 



1 



r/2(k„ 



^ Ceff + 7eff 



k. 



pn 



^/nuj{l - Ccflf) 



Aokpnr]{K) 



'noj{l + C 



cfr; 



n 



2nf 



2nQ 



LU 1 + Ccff 1 + C 

For kz <^ k, we can make the approximation 



(B19) 
(B20) 

(B21) 
(B22) 



tanh"^ A / 1 — 



In 



1 + ^ l-^ 



tanh"^ 



1 + 7' 



(B23) 



1 + 72 ^ V 1 + 7' 2fc2 
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And since we are comparing our result with [8] , we make the following approximation on 
the Coulomb-corrected angle 0^: as the corrections Ak^ and Aky are generally small, we 
can expand to first order in these deviations to write 0^ as a sum of the SFA velocity phase 
and a small correction 6 defined as: 

e tan 6-,, 



tan 5 



where e 



Vx 



1 + (1 + e) tan0^, 
This way we can split the exponential e '^'^^^i'^^i^" )] ; 



(B24) 



(B25) 



A further expansion of (p^ (^^s^"*) "^^^ t)e achieved around At'j^'^ to get: 



exp 



kp - Aoe=°"'""' 



exp 



coAtf^ 










ujAt'P / 


7' I 



^cff - 7eff 



1 + C 



off 



(B26) 



As the probability of escape of the electron in the direction perpendicular to the field is 
exponentially supressed, we can make the approximation k^ ^ kn which gives us: 



kp - Aoe™^"^"' 



-C-(i-C)& + 



= (-1)1™! (^14 
to first order in kz and e{kz 



\m\ 




1+7" 



1 + 



1+7^ 



1 + 



Y J (1 - C2)|m| 



y^-Csgn(m) 



2 

2|m| 



(B27) 



The second term in Eq. (B26) when expanded in powers of kz has a fourth-order depen- 
dence on k/. 



Ceff - 7 



2 

eff 



1 + C 



eff 



C-7^ 

1 + C 



^2(1+72)2 



(B28) 



Finally, we are left with 

2 K2i + 1 {i- \m\)\ 



Wn{S,Uj) = \C, 



Ki\ 



n 47r (£+|m|)! 



(l-(-l)V"")' 



1 + 7' 



2\m\ 



C sgn(m) 



1 + 



-2m 



r 



1+C 



:i-c2)H 

' ,2Wci+2Wc2 



.^^(tanh-y|±J-y|±^ / l+y 2m^K^y>)] (B29) 



d/c, e 



fen 



e + 7' 
kz 



PI' 



±iK 



34 



upto second order in kz- The Coulomb correction is taken out of the integral, on account of 
its extremely weak dependence on the fc^-component of the momentum. The above result is 
valid for all values of i and m. An m-dependent correction due to the Coulomb potential is 
also seen to manifest through its effect on the SFA saddle point t'j'^\ 

To compare with we conisder the case of £ = 1, m = ±1, for which we have P™ (^) = 
— \Jl + ^. To first approximation, we ignore the ^ term in the prefactor, and note that 
since n ^ 1, we can approximate the integral as: 



dKe V I dhe y ^+-^ = K^l- ( ) (B30) 



n f 1+Y 



which gives 



(B31) 



The main difference from Eq. (19) of [8] is incorporation of Coulomb correction, starting 
from the tunnelling region and into the continuum until the electron is registered at the 
detector, and an orbital dependent Coulomb-correction, a result that wasn't expected. 



Eq. (B31) is equivalent to Eq. (70) obtained within time-domain approach. However, 



here we have a result that is valid beyond the optimal momentum, whereas in Eq. ( 70 ) we 



have effectively derived the total ionisation rate summed over all photon orders, which is to 



be compared with Eq. (6) of [8]. For Eq. (B31), further discussion on its range requires a 



knowledge of Ap over all p, and will be considered elsewhere. 
Appendix C: Subcycle ionisation amplitude 

We now consider the case of subcycle ionisation amplitudes in time domain, to replace 
T — 7- 1. The sub-cycle ionisation amplitude is defined as: 

a^{t)=t I cir(p + A(t)|r)^(r,t) (CI) 

J a 

Back-propagating the solution tp{r,T), we can write ip{r,t) as: 
^(r,t)= [ dr' G{r,t]r',T)^{r',T) - I [ dt' [ dr' G{r,t;r' ,t')6{r' - a)B{a,d' ,t') (C2) 

Ja Jt Ja 
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The second term represents that part of the wavefunction that remains bounded within 
the confines of the Coulomb potential near the atom after ionisation. But the wavefunction 
content in that region after ionisation is negligible compared to the current fiux in continuum, 



thus making the contribution from the former almost zero. So we can write equation (CI) 
as: 



apit)=i j rfr(p + A(t)|r) j rfr' G^v^(r, t; r', r)V^(r', T) 

-i(p+A(t))T gj(k+A(t))T 

(27r)3/2 (2^ 



-i(p+A(t))-r i(k+A(t))-r-jk-r' ,t , , 
i dr dr' dk — — g-i/T"!rt^(rL(r;r,p,i)) 



'a J Sl 



Before we can perform the integration on r, we need to address the (r, k)-dependence of 



the Coulomb correction in the above equation. Similar to section |III[ we will expand the 
Coulomb phase term Gc(k, T; r, t) = J^drU (r + dCv]i{C)), about the appropriate sad- 
dle point upto quadratic terms in deviation (a — r^). We need the saddle point for the 
phase term: 

^s^^(r, k, t) = (k - p) ■ r - i dr vUr) (C3) 



Therefore, 



and 



-<SFA 



Vk5^^^ = r- / rfrvk(r) =0 (C4) 



Vr^^^^ = k - p = (C5) 
So the classical trajectory can be written as: 

'c?rvp(r) (C6) 
After expanding the Coulomb phase term Gc{^,T;r,t) about the saddle points (rs,ks) as 



t' 



in section III, we can write the subcycle transition amplitude as: 



Op 



(27r)'^ J J 



Note the subscript p in Gc- the phase term is evaluated for the asymptotic momentum p 
and hence the corresponding momentum shift from this Taylor-expansion Ap = VGc is 
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also evaluated for the asymptotic momentum p and not for the intermediate momentum k 
on which we have to perform the integration. 

We also have achieved the conventional ad-hoc method used to calculate momentum shifts 
[23] . into a rigorous anaylsis arising from stationary points of the EVA action. Following 
our analysis, we first propagate the electron till the detector after ionisation, and to find the 
momentum shifts at any point of time during this motion, we propagate it back through the 
EVA Green's function and thus have information on sub-cycle momentum shifts also. 

We can now write, 

fdTu(r+ rdCMo) = f drul r dCMC)] (C8) 

JT \ Jt J r=r,,k,=n JT \Jt' J 



And we can combine this with: 



dr\j[ l^^JC^m (C9) 



in ap(T) Eq. (36), to get: 




which solves the Coulomb correction for ap{t). The integral on k this time gives k = 
p + Ap(t,T). The r-integral now is straightforward to calculate, and we finally get: 

(cii) 

where we have ignored corrections of the order of OIG^) and greater, which would arise from 
the Coulomb phase and the Coulomb-shifted velocity phase (pl after taking k = p + Ap(t, T). 
Expanding J^*(o) dr Vp^^^^T) upto first order in Ap, it will cancel the spurious term ■ Ap. 
Also p'^ is defined as p'^ = p — Ap{tT\t). The final expression for sub-cycle transition 
amplitude is: 

(C12) 
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Appendix D 



We derive here the result 



pJi{pd(kp,pp)) 



We start from the integral 

Ip^ d(t) dp' p'e'^^""^"^"' (D2) 
Jo Jo 

This integral can be written as: 

I — I del) I dp' p'e^^'^''''' 
Jo Jo 

^ ^ ni=— oo n2=— oo 

= 27r r dp' p'Jo{p'd{kp,Pp)) 
Jo 

In going from step two to three, we first perform the integral over (p and then use the Graf 
generalisation of Neumann summation. The integral over p' is simple 

r*VV„(k,,pJ = ^i!M%f)i (D3) 

Jo "l*^P)Ppj 



Therefore 



Now, by definition 

"277 /"OO 



/'ZTT poo 

{2nfd{kp-pp)= d(j) dp' p's'^^^-P^^p' 
Jo Jo 

p2ii pp 

= lim / #fc / dp' p'e'^^''-^''>P' 
p-^°° Jo Jo 



And we get 



lim 2^ My(k„P,)) 



P-^oc d{kp, pp) 



which is the desired result. 
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Appendix E 



We establish the relation 



kp - Aoe 



where = tan" 



-1 %W 

Vx(t) ' 



We can write 



(pvit) = tan 



tan 



kp sin (f)k — Aq sin out 
kp cos (pk — Aq cos out 



kpc'^'' - Aqc^^ 



i4>k 



Aoe'^^) 



Taking $ = 



pe-^'I'k-Aoe-i'^^ 

(j)v{t) — tan~^ 



i {kpC'-'t'k - Aoe^'^^ + kpe'^'t'k - AqC-^'^*) 
) , we get 



.1 - $ 
1 + $ 



itanh ^ 



1 - $ 
1 + $ 



= - In 
2 



1+1^ 



1-$ 
1+4" 



-- ln$ 
2 



■iujt 



- In , 

2 V^-pe-'"i*fc - AqC- 



Therefore, 



exp 



In 



lOl, 



kpt'^'' - Aqc''^^ 



iojt 



h 1 ,-/(<■',, --■t)\ ™/2 

A;. - A^e^i't'k-'^*) 



(El) 



(E2) 



(E3) 



which is the required result. 
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